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Abstract 

Three methods of probabilistic uncertainty 
propagation and quantification (the method o! 
moments, Monte Carlo simulation, and a 
nongradient simulation search method) are applied 
to an aircraft analysis and conceptual design 
program to demonstrate design under uncertainty. 
The chosen example problems appear to have 
discontinuous design spaces and thus these 
examples pose difficulties for many popular 
methods of uncertainty propagation and 
quantification. However, specific implementation 
features of the first and third methods chosen for 
use in this study enable successful propagation of 
small uncertainties through the program. Input 
uncertainties in two configuration design variables 
are considered. Uncertainties in aircraft weight are 
computed. The effects of specifying required 
levels of constraint satisfaction with specified 
levels of input uncertainty are also demonstrated. 
The results show, as expected, that the designs 
under uncertainty are typically heavier and more 
conservative than those in which no input 
uncertainties exist. 

Introduction 

The aerospace vehicle design process is 
inherently a multidisciplinary design optimization 
(MDO) problem 1 ' 42 . Within recent years, such 
MDO problems have received a growing amount 
of attention from both the engineering and 
optimization communities using both gradient- 
based and nongradient optimization methods. 
Indeed, a quick survey of recent conference 
proceedings and of internet sources reveals 
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literally hundreds of papers addressing various 
aspects of just gradient-based optimization for 
aerospace configurations and its components 

A mission analysis module or discipline is 

usually at the core of an aerospace vehicle 

optimization problem. For example, the major 

contributing disciplines of an aircraft design 

problem (most notably, aerodynamics, structures, 

' * 4 (> 8 C 18 

and performance) may be interacted 

41 42 by using an aircraft mission analysis module. 

Also, constraints formulated within a mission 

analysis module may be used to account for other 

contributing disciplines (such as aircraft layout) 

that are difficult to implement, and for features of 

the vehicle (such as empennage) or the mission 

(such as takeoff and landing) that are not the 

primary focus of the particular design study 41 4_ . 

A particular mission analysis implementation 

known as the Flight Optimization System 

(FLOPS) for aircraft 43 is chosen here as the basis 

for further study. 

In many published optimization studies, the 
inputs to the disciplinary analyses and to the 
multidisciplinary optimization are assumed to be 
precisely known for a given problem; these studies 
are henceforth referred to as deterministic analyses 
and optimizations. In the last few' years interest 
has grown, particularly within the structures 
discipline, in solving problems for which the 
inputs are uncertain 44 " 78 : however, such 
nondeterministic optimizations are relatively 
uncommon in the aerodynamics and related 
disciplines. 

Uncertainties are a prominent aspect early in 
the design process of a new aerospace vehicle, and 
these uncertainties should be accounted for in a 
formal way. The uncertainty in inputs lor 
nondeterministic studies may be due to accepted 
approximations, unmodeled physics, a lack of 
knowledge 79, 8() about some aspect of the problem, 
or errors, such as a lack of precision or 
repeatability in measurement, or blunders 
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attributed to the user or their processes. In these 
cases, the uncertain inputs may be considered 
random variables that take on some prescribed 
distribution of values rather than a single precise 
value. The uncertain input value distribution may 
be one of numerous popular distributions 
categorized by various statistical and probabilistic 
resources* 0, 81 , or it may he some uncategorized, 
perhaps even unknown, distribution that can only 
be approximated. These studies, in which input 
variables are assumed random and are drawn from 
a prescribed distribution of values, are henceforth 
referred to as nondeterminislic analyses and 
optimizations. The potential for input uncertainty, 
expressed to disciplinary analyses via distributions 
rather than single values, introduces a level of 
uncertainty in the resulting analysis outputs and 
raises the need to consider output distributions. 

A wide variety of probabilistic uncertainty 
propagation techniques exist 82-85 . The uncertainty 
propagation techniques generally fall into one of 
six categories: simulation methods, importance 
sampling techniques, first-order reliability 
methods, second-order reliability methods, 
response surface methods, and method of moments 
techniques. This paper considers the method of 
moments, a Monte Carlo simulation technique, and 
a nongradient simulation search method. The 
Monte Carlo simulation technique is only used for 
comparison with the method of moments 
approximation. 

The application of probabilistic methods requires 
the definition of ( 1 ) one or more random input 
variable probability models or distribution types, 
(2) one or more response models that describe the 
physics, process, or rules which govern the system 
behavior, and (3) one or more models that predict 
the outcome of an event; these predictive models 
are generally called limit states. The intent of each 
of these uncertainty propagation techniques is to 
evaluate a multidimensional probability integral 
over a multidimensional surface known as the 
limit state. However, in practice, evaluating this 
multidimensional integral in closed form is 
problematic for several reasons: ( I ) the joint pdf is 
generally not known, (2) the boundary over which 
the integral is to be evaluated (the limit state) is 
generally not known, and (3) even when the pdf 
and limit state are known, the multidimensional 
integral itself is difficult to evaluate. As a result, 
various uncertainty methods, each with differing 
computational features, and levels of accuracy and 
efficiency, have been proposed and developed to 
circumvent these difficulties. 

The random variable probability models for 
this study are chosen for convenience from among 


many possible distribution types; for example, the 
variable probability models could be described by 
normal, lognormal, Weibull, uniform, or beta 
distributions. For each distribution type, the pdf 
describes the probability that a certain value of the 
random variable will occur, plotted as a function 
of the range of possible values that can be assigned 
to the random variable. The shape of these 
distributions is generally described analytically by 
at least two parameters, including the mean value 
(denoted herein by variables with an overbar) and 
the standard deviation (a), which is a measure of 
the dispersion of the random values about the 
mean value. Some distributions may require more 
than two parameters to be described, but the mean 
value and standard deviation are sufficient for the 
distributions considered herein. The standard 
deviation is the product of a random variable mean 
value and the more commonly chosen coefficient 
of variation (c.o.v.). The reader should note that 
the normal, lognormal, and Weibull distributions 
are unbounded in at least one direction, whereas 
the uniform and beta distributions are bounded in 
both directions. In this paper, only normal 
distributions are used for input variables. 

Each pdf has an associated cumulative 
distribution function (cdf)describing the 
probability that the value of a random variable is 
less than or equal to some prescribed value taken 
from the total range of possible values. Two or 
more random variables may be correlated 
(dependent, or unrelated but changing together), or 
uncorrelated (totally independent); in this paper, 
only truly independent random variables with 
normal distributions are considered. 

A given uncertainty analysis might yield one 
or more of the following results: identification of a 
single most probable point (mpp), or a locus of 
mpp, at which a certain event might occur; 
computation of the reliability index and its 
gradients with respect to the mean values and 
standard deviations of the random input variables 
at the MPP; computation of the parameters 
describing a random output pdf or cdf given one or 
more known random input pdf and cdf; 
approximation of the output pdf or cdf shape 
without the formal calculation of the parameters 
associated with such a distribution; or the 
accommodation of uncertainty in random input 
variables without approximating the pdf and cdf. 
The method of moments used in this paper 
accommodates uncertainty in random input 
variables and produces estimates for the mean 
value and the standard deviation of the output 
variables. These estimates are accurate only if the 
resulting output distribution is normal. The Monte 
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Carlo simulation method can approximate the pdf 
of an output function for any input distribution, 
without actually computing the parameters needed 
to describe the pdf shape, but requires many 
function evaluations for accuracy. In contrast, the 
nongradient simulation search method (ssm) can 
efficiently produce very accurate representations 
of the pdf and cdf for any distribution types, as 
well as identifying the mpp, the reliability index, 
and the gradients of the reliability index w'ith 
respect to the random input variable mean value 
and standard deviations at the mpp. 

An uncertainty propagation and 
quantification demonstration for a one- 
dimensional inviscid aerodynamic example 
problem, using normal distributions for two input 
random variables, was developed in a recent paper 
by Putko et al . 86 Parameters for the output 
distributions were computed by the method of 
moments and verified by Monte Carlo simulations. 
Other work in aerodynamic uncertainty 
propagation and quantification has since been 
published involving the extension to two- 
dimensional airfoil design under uncertainty 
and the propagation of uncertainty through 
aerodynamic viscous and turbulent flows via the 
application of random field theory * 0 >(! . The w ? ork 
of Putko 86 was also the basis for similar 
uncertainty propagation work in a three- 
dimensional aerodynamic/slructural interaction by 
Gumbert et al. 91 . The Putko and Gumbert work 
closely followed the work of others from the 
structures discipline, most notably 92 93 , which 
showed that a statistical first-order second-moment 
(FOSM) method and automatic differentiation 
(AD) could be used to efficiently propagate the 
input uncertainties through finite element analyses 
to approximate the output uncertainty. An 
integrated strategy for mitigating the effect of 
uncertainty in simulation-based design is presented 
in 93 ; it consists of uncertainty quantification, 
uncertainty propagation, and robust design tasks or 
modules. Two approaches are developed in 9 to 
propagate uncertainty through sequential analysis 
codes: an extreme condition approach and a 
statistical approach. The latter approach can be 
efficiently implemented using the FOSM 
approximation and sensitivity derivatives (SD) as 
was done by Putko. 

The response model for this demonstration is 
the FLOPS mission analysis code, which includes 
equations to compute the aircraft weight, range, 
cost, noise metrics, and other readily available 
performance constraints. The method of moments 
was implemented in the FLOPS aircraft mission 
analysis program for several classes of statistically 


independent, normally-distributed, random input 
variables and several classes of random output 
variables, noted above. However, due to concerns 
about the accuracy and applicability of the method 
of moments for this code, only uncertainty 
propagation results for examples with two 
uncertain input variables and one uncertain output 
variable are shown in this paper, Diflerent levels 
of input uncertainty and required constraint 
satisfaction are imposed. The effect of uncertainty 
on the design point, compared with a deterministic 
design, is noted. Output distributions from the 
deterministic code are compared with Monte Carlo 
simulations. Sample results from the ssm 
technique of the UNIPASS tool are also 
shown. 

The problem (including the FLOPS mission 
analysis code, the aircraft, its mission, design 
variables, objective, and constraints) was chosen 
for this uncertainty demonstration because this 
particular mission analysis code executes quickly, 
which enables some level of validation with Monte 
Carlo techniques. Furthermore, the FLOPS code 
has been shown to be very amenable to 
processing 42 by the Adifor automatic 
differentiation tool IUOOO:> , which enables very 
efficient computation for thousands of derivatives. 
Unfortunately, the FLOPS code also produces a 
large number of failed analyses (illustrated 
subsequently) when executed for a series of related 
cases, as might be done during Monte Carlo 
simulation. The frequent failures of the analysis 
code result in a design space that, to an 
optimization or uncertainty propagation tool, 
appears to be discontinuous. This situation is 
actually a common feature of many complex 
analysis codes. It is also a feature that is not well- 
handled by many commercially developed tools. 
But the discontinuous nature of the FLOPS design 
space serves as a good example to illustrate the 
benefits of a nongradient technique available in the 
UNIPASS 1 M tool. In this paper, all work related to 
and using the UNIPASS™ tool was performed by 
the PredictionProbe, Inc. experts, whereas the 
method of moments work was performed by the 
NASA civil servant. 


^ The use of trademarks or names of manufacturers 
in this report is for accurate reporting and does not 
constitute an official endorsement, either 
expressed or implied, of such products or 
manufacturers by the National Aeronautics and 
Space Administration. 
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Approach and Methods 

The FLOPS Mission Analysis Code 

The FLOPS code 43 is a multidisciplinary 
system of computer programs for conceptual and 
preliminary design and evaluation of advanced 
aircraft concepts. It consists of nine primary 
modules: weights, aerodynamics, engine cycle 

analysis, propulsion data scaling and interpolation, 
mission performance, takeoff and landing, noise 
footprint, cost analysis, and program control. 

The FLOPS code may he used to analyze a 
point design, to parametrically vary certain design 
variables, or to optimize a configuration with 
respect to numerous design variables using 
nonlinear programming techniques. A variety of 
configuration, mission performance, noise 
abatement, and engine design variables and a 
composite objective function (including the 
minimum gross weight, minimum fuel burned, 
maximum range, minimum cost, and/or minimum 
NO* emissions) are provided to allow for 
simultaneous optimization of the aircraft 
configuration, engine cycle, and size. 

Two example cases, distributed with the 
FLOPS code, will be used for the uncertainty 
propagation demonstrations in this paper. The first 
is a subsonic aircraft transport design case. The 
second is a supersonic aircraft transport design case. 

The AD1FOR Automatic Differentiation Tool 

The ADIFOR 1(KM0 * software package is a 
tool for the AD of standard FORTRAN 77 
programs. Given a FORTRAN 77 source code 
and user-specified dependent and independent 
variables, Adifor will formulate exact derivatives 
(via repeated and systematic application of the 
chain rule of calculus) and generate new f 
FORTRAN 77 code. The new code includes 
original function evaluation, augmented with code 
that computes the partial derivatives (gradient) of 
the specified dependent variables with respect to 
the specified independent variables. 

The ADIFOR 2.0 software package 102 
provides a production-quality AD environment 
that can compute derivatives by the forward 
(direct) mode of AD. In the forward mode of AD, 
the gradient code execution time and memory are 
usually proportional to the number of independent 
variables: this technology is best suited to 
problems in which the number of dependent 
variables is greater than the number of 
independent variables. 

The ADIFOR 3.0 software package 10 ’ 105 
includes both forward and reverse (adjoint) modes 
of AD for first derivatives, and three forward 


mode options for computing second derivatives. In 
the reverse mode of AD, the gradient code 
execution time and memory are usually 
proportional to the number of dependent variables; 
adjoint technology is best suited to problems in 
which the number of independent variables is 
greater than the number of dependent variables. 

Previous studies with the FLOPS code 42 
indicate that this program requires only a few 
minor changes to correct nonstandard FORTRAN 
77 coding to enable ADIFOR 2.0 processing of the 
code to compute first derivatives. The resulting 
gradient code was also found to be exceptionally 
efficient in computing thousands of first 
derivatives via the forward mode of 
differentiation; it was almost as quick as the 
original code execution. Both ADIFOR 2.0 and 
ADIFOR 3.0 were applied to the FLOPS code 
during the course of this study; in fact, several 
"bugs" in the ADIFOR 3.0 package and in the 
FLOPS code were identified and corrected as a 
result of this work. However, only results using 
ADIFOR 2.0 are shown in this paper, 

ADIFOR Application to FLOPS 

The independent variables for differentiation 
were selected from among the possible 
configuration, mission performance, and noise 
abatement design variables input to FLOPS. The 
independent variables for differentiation also 
included representative elements in the 
aerodynamic and propulsion data provided to the 
code from external sources. The dependent 
variables were selected from numerous outputs 
from the FLOPS analysis, including the composite 
design objective, the takeoff gross weight, the 
vehicle life cycle cost per unit, a noise metric, an 
emissions metric, and seven typical constraints 
from a menu of nineteen that could be activated for 
the problem. 

Since the FLOPS analysis is embedded 
within the optimization, extra care was taken to 
ensure that correct derivatives were obtained for 
the analysis portion of the code during both 
analysis and optimization modes of the code 
execution. This entails allowing for differentiation 
with respect to the input design variable values 
(used for analysis mode, and as the starting point 
for optimization) and with respect to the local 
design variable values used in the analysis module 
when embedded within the optimization mode. It 
also required special handling of the derivative 
activation ("buddy variable") sites within the code, 
the derivative seeding to an identity matrix, and 
the zeroing of certain iterated variables to ensure 
that derivative objects were not contaminated 
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during analysis N with information from analysis 
N- 1 . It is expected that any analysis embedded 
within an optimization will require similar 
techniques to ensure proper derivative evaluation. 

The computed derivatives were only used to 
provide gradients needed for the uncertainty 
augmentation, described subsequently. ADIFOR 
could also be used to provide the derivatives 
needed for the gradient-based optimization within 
FLOPS; however, that was not done here. Higher 
order derivatives could also be obtained via 
ADIFOR differentiation to enable the formulation 
of higher order uncertainty corrections (more 
accuracy) and higher moment approximations 
(more information about the output distributions). 

TM 

The UNIPASS Software 

The UNIPASS™ software, developed by 
PredictionProbe, Inc., is a general-purpose Unified 
Probabilistic Assessment Software System that 
performs complex probabilistic analyses. The 
UNIPASS™ software can be utilized 
independently as a stand-alone software engine, 
and/or integrated with deterministic software tools. 
UNIPASS™ provides the basis for modeling 
uncertainties, constructing probabilistic predictive 
models, computing probabilities, identifying most 
likely outcomes, providing sensitivity data, 
identifying key drivers, analyzing risk, and 
performing sensitivity analysis; deterministic 
software tools may be integrated to provide the 
computational framework for constructing 
complex deterministic process models. The latest 
version. UNIPASS™ 4.2, offers an advanced 
graphical Windows environment. 2-D and 3-D 
graphic functions, four problem types (component, 
serial, parallel, and general), three analysis types 
(probability, inverse probability, and cdl/pdt 
analysis), six categories of probabilistic methods, 
thirty-seven probability distribution types that can 
be used to define any type of random variable, four 
classes of random variables, three ways to 
interface with any commercial and/or in-house 
software tools, eleven gradient-based mpp 
identification methods, and one nongradient 
simulation-based search method (SSM) that finds 
the mpp for discontinuous, nondifferentiable limit 
state functions. In general, the mpp represents the 
most likely values of the random variables at 
which the critical or significant condition of the 
user-defined event will occur. In engineering, a 
critical condition may be an undesirable event 
such as component failure or instability, or a 
desirable event such as extended component life or 
mission success. Some software products use only 


a gradient-based algorithm to identify an mpp. 
However, those algorithms are limited to 
continuous and differentiable variables, and cannot 
handle the more common engineering tasks, such 
as the FLOPS examples studied herein, which 
involve discontinuous limit state functions. 

Uncertainty Propagation 
Method of Moments 

Uncertainty propagation is accomplished by 
using various orders of approximations to the 
various statistical moments; this is a logical 
naming convention for the uncertainty propagation 
technique results from a given choice of the order 
of approximation and the statistical moment to be 
used. Only the FOSM approximation is used lor 
results in this paper. 

In this study, the effects of uncertainty in two 
input aircraft design variables are considered tor 
the purposes of illustrating uncertainty propagation 
through the FLOPS code and design under 
uncertainty. For the present demonstration (and 
following the derivation in 8Sand89 ), these input 
variables are assumed statistically independent, 
random, and normally distributed about a mean 
value. These assumptions simplify the 
implementation and help quantify the input 
uncertainties. The assumption of the variables 
being statistically independent is not required; 
correlation between the variables can be easily 
accounted for within the formulation at the cost of 
more computational work. For non-normal input 
distributions, the method of moments corrections 
are only approximately correct. In this case, the 
increase in accuracy, gained from considering 
additional terms in the Taylor series expansion, 
could very well be offset by the approximation 
error due to the non-normal nature of the input 
variable. 

Givena vector B = {b r ...,b„} with n 
independent input random variables, /? (for i 
= /,/?>, mean values B = j . standard 

deviations a b = {o bj c bn ( , and random output 

function F, first-order (FO) and second -order (SO) 
Taylor series approximations to the function are 
given in generic form by 

FO: 

F(B) = F(B) + ±¥-{b,-b l ) (l: 
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SO: 




2\~^dbdb : 1 ' ' H ' '' 


( 2 ) 


where the first and second derivatives are 
evaluated at the mean values,/?, and B . 

One must then obtain the expected values for 
the mean (first-moment) and variance (second- 
moment) of the output function, F, which depend 
on the derivatives of Fwith respect to the uncertain 
input variables and input variances, o b \ The 
expected values of a random function are obtained 
from the integration of the product of the function 
itself and the imposed pdf. For normally 
distributed input values, the pdf is symmetric 
about the mean value. Thus, the expected value of 
an odd function with a normal input distribution 
involves the integration of the product of an even 
and odd function, which is zero. Likewise, the 
expected value of an even function involves the 
integration of the product of two even functions, 
which is nonzero. The mean value of the output 

function F and standard deviation o F , are 
approximated (as in 8<> ) as 


FO: 

F = F[h) 


I 


dF 

db, ^ 


(3) 


particular sample cases distributed with the 
FLOPS code. The first example uses the inputs for 
the FLOPS five -design variable subsonic transport 
design (xfp2.in). The input file is modified to 
allow only the variables THRUST (the maximum 
rated thrust per engine, in pounds force), and SW 
(the wing reference area, in square feel), to be 
active design variables; upper and lower bounds 
have also been specified for these design variables. 
The optimization objective is specified to be the 
aircraft gross takeoff weight. The seven possible 
aircraft performance constraints, normally 
activated with this sample problem, are used. 

These include the aircraft required range (which is 
held fixed for this problem), the approach speed, 
the takeoff and landing distances, the approximate 
missed-approach and second-segment climb 
gradients, and the excess fuel. The Broyden- 
Fleteher-Goldfarb-Shano (BFGS) optimization 
method (the default among several optimization 
methods available within the FLOPS code) was 
used to solve this problem. In the FLOPS 
implementation of this optimization method, a 
composite objective function is minimized. The 
composite objective function is composed of the 
true objective augmented with a highly nonlinear 
penalty function that grows rapidly as the design 
variables approach their upper or lower bounds, 
and as constraints become active. 

Robust Optimization 

The form of the uncertain objective follows 
from the development in 86 , with an adaptation to 
the current optimization problem. The FOSM 
expression for the current uncertain objective is 


SO: 


L r { l>\ ^ V ^ ^ 

F - F + — > — — a; 

1 ; 2 \ttdb; ■ 




dF ) 

+JLy y 

' d : F ) 

k '■> 

■yi 4A 

- ■ ;=l ; = l 

dbdb h l> 

V ’ / 


: (4) 


where the first and second derivatives are again 
evaluated at the mean values, b, and B . Note in 
Eq. (4) that the second-order mean output F is not 
at the mean value of input/? ; a shift in the mean 
value of the uncertain output function occurs due 

to the specified input uncertainly, i.e., F * F^B j . 


Deterministic Optimization 

For simplicity, a demonstration of 
deterministic optimization is derived from two 


Obj = 


. Obf 

i.o+ a ° hl } 

V 

Obf ) 


(5) 


where 




dO hi 

"d/P 


( 6 ) 


and where the o b are the known standard 

deviations for each of the random input variables, 
as noted previously. The adaptation of the 
uncertain objective provides a composite objective 
function, similar in magnitude and functionality to 
the original objective function of the FLOPS code 
but augmented with uncertainty effects. Since the 
uncertainty correction to the objective function is 
small compared with the uncertainty correction to 
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the constraints, it may be possible to neglect the 
objective uncertainty correction altogether; this is 
proposed as the subject of further research. 

Again, following the development in S( \ the 
forms of uncertain constraints are written as: 

e = e +k(J . m = \J (7) 

c> m o m x,„ 


where 



m = 1,7 


( 8 ) 


and where k is the required probabilistic constraint 
satisfaction in units of input variable standard 
deviations. For this example, values of k ranging 
from zero (constraint satisfied with 50% 
probability for a normal distribution) to three 
(constraint satisfied with 99.9% probability for a 
normal distribution) were considered. Values of 
the random variable <7 h were each computed as 

the product of a c.o.v. and the random input 
variable mean value, both input to FLOPS 
program, for each variable. The c.o.v. values were 
chosen to be the same for all the random input 
variables in a given problem. 

It is surprising that tiny numerical differences 
between the deterministic and robust solutions 
(effects of uncertainty) could be observed for 
extremely small values of the input c.o.v. (i.e., for 
c.o.v. of order 10 20 ). These solutions are 
“deterministic" from a practical point of view, but 
exhibit numerical behavior that can only be 
attributed to the imposed uncertainty corrections. 
However, for "small enough" values of c.o.v. (i.e., 
for c.o.v. of order l O' 30 ), no differences from the 
purely deterministic solution could be discerned, 
as expected. 

The method of moments formulation 
described above requires derivatives to be 
evaluated. If these derivatives were to be 
computed by finite difference approximations, it 
would be necessary to find two or more successful 
function evaluations near the point of interest for 
evaluating the derivative. For the example 
problems chosen, this requirement was sometimes 
difficult to meet because of frequent failures of the 
analysis module, which led to discontinuous 
objective and limit state functions, although the 
design physical space of interest was locally (or 
piecewise) continuous. The discontinuous nature 
of the design space with the FLOPS code also 
caused problems for the UNIPASS 1 M tool, 


necessitating the use of the nongradient SSM. 
However, the use of automatic differentiation to 
compute the derivatives made the method of 
moments formulation more successful than might 
be expected for examples with numerical 
discontinuities because the derivatives were 
evaluated analytically at each point of interest via 
the chain rule. Thus, finding successful 
neighboring points to a successful analysis point 
was not an issue with the current method of 
moments uncertainty formulation. 

Two restrictions that are more serious do 
arise from using the method of moments: ( 1 ) the 
formulation is only valid for normal input and 
output distributions, and (2) the function gradients 
are evaluated at the mean value point, rather than 
at the mpp. Thus, the method of moments 
computations are generally not very accurate away 
from the mean value of the random function and 
are even less accurate for non- normal 
distributions. In particular, the method of 
moments is expected to be significantly inaccurate 
for predicting very low probability of failure 
points (e.g., in the tail region of a normal 
distribution) and for highly skewed output 
distributions. The use of AD to evaluate gradients 
offers two possibilities to address this issue: ( I ) if 
higher order derivatives are available, better 
approximations can be achieved for low 
probability failure points by constructing Taylor 
series approximations that better represent those 
evaluation points far from the mean value, and (2) 
if the mpp can be identified by other means, the 
functional derivatives can be obtained at the mpp 
via AD without the need to find successful points 
neighboring the mpp. 


Simulation Search Method (SSM) of 
UNIPASS™ 

The FLOPS optimization routines, like many 
other optimization algorithms, may not produce 
acceptable results for a given set of input data, and 
therefore will result in a highly discontinuous 
response surface. When probabilistic analysis is 
performed, such discontinuity in the response 
surface cannot be handled with standard first- and 
second-order reliability methods or other 
approximate and more efficient probabilistic 
techniques. Furthermore, since the desired level of 
reliability is usually greater than 0.99999, the 
application of Monte Carlo simulation becomes 
impractical. Therefore, an efficient first- and 
second-order reliability method that could solve 
such nondifferentiable discontinuous problems. 
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such as the one provided by UNIPASS™, 
becomes highly desirable. 

The UNIPASS™ nongradient SSM was 
applied to the FLOPS code for the sample problem 
based on the subsonic transport design problem, 
again modified to allow only the THRUST and 
SW to be active design variables. Application of 
the SSM also allowed for: ( I ) determination of the 
aircraft weight PDF for a probability range from 
0.01 to 0.99, (2) determination of the aircraft 
weight CFD and the maximum aircraft weight due 
to input uncertainties to a probability of 0.9999, 

(3) mpp identification, and (4) computation of the 
gradients of the reliability index. 

Results 

Method of Moments 

Starting from the subsonic transport 
optimization problem (“xflp2.in") test case 
distributed with the FLOPS code, the input file 
was first modified so that only THRUST 
(maximum rated thrust per engine) and SW (the 
reference wing area) in Namelist SCONFIN were 
active design variables. These variables are the 
primary aircraft sizing variables. The modified 
input file was then used for both deterministic and 
robust optimizations. The robust optimizations had 
various levels of input uncertainty for the active 
design variables and various levels of required 
constraint satisfaction, both specified in auxiliary 
input file to FLOPS. For the first example case, 
the specified uncertainty corresponds to a c.o.v. of 
5 % for each of the two input variables. Since the 
mean value of THRUST for this problem is 47,500 
lb and the mean value of SW for this problem is 
2272 ft 2 , one standard deviation (a) is 2375 lb, and 
1 1 3.6 fr for the two variables, respectively. 

Figure 1 is a simplified aircraft sizing contour 
("thumbprint") plot illustrating the design space 
near the deterministic and robust optimizations for 
this case. The x-axis of the figure shows a 
normalized value of the maximum rated thrust per 
engine, ranging from 20,000 to 60,000 pounds 
force; the y-axis of the figure shows the 
normalized value of the reference wing area, 
ranging from 1000 to 5000 ft 2 . 

The figure is simplified from typical 
thumbprint plot in that only the active constraint 
violation boundaries for the specified problem are 
shown, rather than multiple contours representing 
various values of the limit slate function. The 
constraint violation boundaries were interpolated 
from a parametric variation of the two variables 
with nine equally spaced points in each direction. 
Also shown in figure I are the locations of final 


design points for the deterministic and three mean 
robust optimizations. In this example, both 
FLOPS constraint 2 (the upper limit on approach 
speed) and FLOPS constraint 5 (the lower limit on 
missed approach climb gradient) were active for 
the deterministic optimization. As expected, the 
deterministic design point is found at the 
intersection of the two active constraints. The 
optimization path was almost entirely w ithin the 
feasible region in the figure. 

The three robust optimization points (labeled 
A - 1 , k = 2, and k = 3 in the figure) correspond to 
imposed constraint satisfaction margins of I, 2, 
and 3 standard deviations about the mean value of 
the deterministic solution. The offset in the robust 
design points from the constraint violation 
boundaries is proportional to both the imposed 
input uncertainty and the gradient of the constraint 
with respect to the uncertain design variables. 
Figure 1 also shows that the robust optimization 
with k = 1 enforces a greater margin of 
satisfaction for both constraints than does the 
deterministic optimization. Similarly, each of the 
robust solutions enforces greater constraint 
satisfaction, with increasing values of A than either 
the deterministic solution or the robust 
optimizations with smaller values of A. For the 
deterministic optimization, the constraint 
satisfaction with respect to single constraint 
violation is only 50% probability for an output 
normal distribution; for A = I, this probability 
increases to about 84%; for A = 2, the probability 
is about 97.7%; and for A = 3, the probability rises 
to about 99.9%. 

Simultaneously (not shown in the plots), the 
aircraft weight increases with increased constraint 
satisfaction from the deterministic value of 
21311 0.5 lb to a value of 2 1 7052.6 lb for A = I , to 
a value of 220222.0 lb for A = 2, and to a value of 
223443. 1 lb. for A = 3. The weight for A = 3 is 
about 5% higher than the deterministic solution. 

Similar results to those shown in figure I are 
presented in figure 2, this time for a c.o.v. of 10%. 
The relative offset from the constraint violation 
boundaries grows in proportion to the increased 
input uncertainty imposed on this problem, relative 
to the previous example. Simultaneously (not 
shown in the plots), the aircraft weight increases 
with increased constraint satisfaction. The 
deterministic solution has a weight of 213] 10.5 lb; 
A = I has a weight of 2202 1 8.9 lb; k-2 has a 
weight of 227079.4 lb; and A = 3 has a weight of 
234092.7 lb (about 10%^ greater than the 
deterministic weight). 

Generally, similar results to those shown in 
figures 1 and 2 are shown for the supersonic 
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transport optimization problem (“xflp3.in”) lest 
ease in figures 3 and 4 with violation boundaries 
constraints 2 and 3 (takeoff distance). In this case, 
the design under the uncertainty problem was 
much more difficult to solve in a robust design 
mode: the design path was almost entirely in the 
infeasible region of the figure, meaning the 
uncertainly correction contributed substantially to 
the nonlinear penalty function. The example 
exhibited much greater sensitivity to smaller levels 
of input uncertainty and to smaller levels of 
variation in local gradients calculated lor 
objectives and constraints at various points in the 
design evolution. The deterministic design point 
appears to be caught in a corner of constraint 2, 
rather than at the intersection of the two sometimes 
active constraints. In many cases, the optimization 
path was significantly different under small levels 
of input uncertainty than for the deterministic 
optimization. Figure 3 shows the results for only 
0.05% variation of the THRUST and SW. The 
same general behavior is observed as in figures 1 
and 2. but figure 3 shows much greater levels of 
sensitivity to the level of input uncertainty. Figure 
4 shows the results for 0.08% variation in the same 
two input variables. 

Results from a 5000 sample Monte Carlo 
simulation, centered at the deterministic design 
point of figures 3 and 4 and based on 106 , are shown 
in figure 5. Of the 5000 requested FLOPS 
analyses. 447 (8.94%) failed to produce an answer. 
Over the course of preparing this paper. FLOPS 
analysis failure rales ranging from 8.2% to 88% 
were observed during various NASA civil servant 
attempts to perform Monte Carlo simulations, 
depending on the various parameters chosen to 
guide the Monte Carlo simulation. The average 
analysis failure rate from seven Monte Carlo 
attempts was 28% , which agrees closely with the 
FLOPS analysis failure rate observed by the 
experts from PredictionProbe during the course of 
their studies. The FLOPS gross weight response 
(labeled Prob in the figure) was converted to a 
standard normal space for comparison with a 
standard normal distribution (labeled Norm in the 
figure). The figure shows that the FLOPS output 
distribution departs significantly from that of a 
normal distribution, suggesting a highly nonlinear 
response from the FLOPS code for this supersonic 
transport design case. The mean value shifts 
significantly from that of a standard normal 
response and even a bimodal response pattern is 
observed, indicating that the lower aircraft weights 
could be found " just around the corner” from the 
deterministic optimization point. This result might 
be expected after examining the location of the 


deterministic design point in figures 3 and 4 at a 
corner of a single constraint in the design space. 

Despite the obviously non-normal output 
behavior of the FLOPS code for this example, as 
depicted in figure 5, on-going research suggests 
that the failure rates with respect to constraint 
violation predicted by the method of moments are 
similar in magnitude to those predicted by Monte 
Carlo simulations centered about the various 
design points depicted in figures 1-4. More 
research into the accuracy of the method of 
moments, perhaps utilizing AD to compute higher 
order approximations to higher moments, might be 
warranted in this case. 

UNIPASS™ Example 

A pdf/cdf analysis was performed for the 
subsonic aircraft design problem using the SSM to 
determine the aircraft weight distribution for the 
probability range from 0.01 to 0.99 (the middle 
98%' of the weight distribution). The SSM allows 
for identification of the mpp and calculation of the 
reliability and sensitivity data for nondifferentiable 
discontinuous problems. Assuming THRUST and 
SW to be normally distributed with mean values of 
34405. 1 1734 and 2054. 19523, respectively, 
figures 6 and 7 depict the pdf and cdf of the 
aircraft weight for the deterministic optimization 
shown in figures 1 and 2, for c.o.v. = 5% and 10 
%, respectively as determined by the UNIPASS IM 
tool. Note that the output distributions are not 
normal distributions, which implies that the results 
shown in figures 1 and 2 tor the method of 
moments uncertainty propagation may also be 
suspected of inaccuracy. 

The SSM is much more accurate than a 
comparable Monte Carlo simulation, for a given 
number of analyses, at predicting the output 
distribution cdf/pdf. For example, more than 
40,000 successful runs would be needed to 
produce an accurate cdf point for a probability 
level of 0.99 with c.o.v. = 0.05% . Furthermore, 
while a 0.99 probability level may be adequate to 
predict the potential overweight, it is far from 
adequate for reliability estimates or risk 
calculations, which often require success 
probability levels of 0.99999+. In such cases, 
40,000,000+ successful runs are needed to predict 
probabilities with c.o.v. = 0.05%, which is 
impractical for real world applications. 
Alternatively, an SSM approach may be used to 
develop the cdf/pdf of the aircraft weight 
significantly fewer runs. In this case, the number 
of the maximum number of runs is determined by 
the SSM based on a predefined tolerance level. 
Using SSM, 457 executions of FLOPS produced 
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263 successful runs (42.4% failure rate), which 
were used to identify an mpp, to produce the cdf 
plots (figures 6 and 7), and to generate sensitivity 
measures (figures 8 and 9). The total number of 
successful runs required to produce accurate 
values for probability values of 0.99999+ is not 
expected to increase significantly. 

Sample sensitivity data from UNIPASS 1 M 
showing the sensitivity of the reliability index with 
respect to the mean values of the input random 
variable, evaluated at the mpp and scaled by the 
random variable mean values, are depicted in 
figures 8 and 9. for c.o.v. = 5% and 10%, 
respectively. Significantly, these figures indicate 
that there are three distinct regions in the design 
space: ( 1 ) negative sensitivities for both SW and 
THRUST, (2) positive sensitivity for THRUST 
and negative sensitivity for SW. and (3) positive 
sensitivities for both variables. Additional studies 
are needed to understand the implications and 
limitations of this sensitivity information for the 
given problem before drawing any conclusions. 

For example, recall that these two variables were 
deemed to be uncorrelated (without discussion), 
when, in fact, greater SW would logically require 
greater THRUST. Subsequent studies should fully 
consider any possible statistical correlation 
between the uncertain variables. It should also be 
emphasized that this probabilistic study was 
performed using only two uncertain variables, 
whereas the potential number of deterministic and 
uncertain variables could be significantly more for 
many cases of interest. 

Conclusions 

The FLOPS aircraft mission analysis and 
optimization code was successfully augmented 
with approximations to the first-order second- 
moment probabilistic uncertainty propagation 
terms for the objective and potentially active 
constraints. Two input variables that substantially 
contribute to aircraft shape and sizing were 
assumed to be uncertain in two separate test cases: 
a subsonic transport design and a supersonic 
transport design. These variables were assumed to 
be statistically independent and to take on random, 
normally distributed input values about a mean 
value. Gradients required for uncertainty 
augmentation were obtained by using automatic 
differentiation applied to the code. 

Results from two deterministic optimizations 
and from several designs under uncertainty, with 
various amounts of imposed uncertainty for two 
input design variables, were presented. For the 
subsonic transport design case, input uncertainties 
of 5% and 10% of the mean value of the uncertain 


input variables were considered. Results were also 
shown for increasing amounts of required 
constraint satisfaction. As expected, the weight of 
the aircraft increases in all cases from its 
deterministic value. The amount of weight 
increase was proportional to both increasing 
amounts of uncertainty and to increasing amounts 
of required constraint satisfaction specified in the 
optimization problem. 

For the subsonic transport design problem, 
the output probability density function 
distributions computed by a commercially 
distributed uncertainty propagation tool are non- 
normal in shape, indicating a nonlinear response 
from the code for which the method of moments is 
known to be inaccurate. For the supersonic 
transport design case, uncertainties of only 0.08% 
of the mean value of the uncertain input variable 
were considered. They produced levels of output 
uncertainty similar to those for the subsonic- 
transport example with 10% input uncertainty. 

The weight distribution also departed significantly 
from that of a standard normal distribution, 
making the application of the method of moments 
highly questionable for this case. Although the 
theoretical accuracy of the uncertainty propagation 
results obtained with the method of moments is 
questionable for the cases shown, in practice, the 
method may still yield reasonable approximations 
for the constraint failure rate. Further research into 
the accuracy of the method of moments 
approximation is recommended for both cases. 

The FLOPS code was also successfully 
integrated and analyzed with a commercially 
distributed probabilistic assessment software 
system to identify aircraft weight probability 
density function, cumulative density function, and 
maximum weight to a high probability level, and 
the most probable point (mpp) of failure with 
respect to the imposed design constraints. Since 
the analysis module under consideration exhibited 
significant numerically discontinuous behavior, 
the mpp was found using a nongradient simulation 
search method. Once the mpp was known, an 
inverse probability technique was applied to 
compute sensitivities of the reliability index with 
respect to the mean value of two random input 
variables. The probability density function and 
cumulative density function of the random output 
variable were also computed. 
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Fig. 1 . Design points and active constraint 
boundaries for deterministic and robust 
optimization, c.o.v. =5%. 
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Fig. 2. Design points and active constraint 
boundaries for deterministic and robust 
optimization, c.o.v. = \0 r /(. 



Fig. 3. Design points and active constraint 
boundaries for deterministic and robust 
optimization, c.o.v. = .05 r /< . 



Fig. 4. Design points and active constraint 
boundaries for deterministic and robust 
optimization, c.o.v. = .08 ( 7c. 



Fig. 5. Aircraft gross weight distribution plot, 
supersonic, two variable, deterministic 
design case. 
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Fig. 6. Aircraft gross weight distribution plot, subsonic, 
two variable, deterministic design case, c.o.v. = 59k 



two variable, deterministic design case, c.o.v. = 10%. 
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Limit State 1 


Fig. 8. Reliability index sensitivities with respect to random input variable mean values, 
scaled by the random input variable mean value, c.o.v = 5 %. 



Limit State 1 


Fig. 9. Reliability index sensitivities with respect to random input variable mean values, 
sealed by the random input variable mean value, c.o.v = 1 ()9r. 
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